Aligning rectilinear images in 3D through projective registration and calibration

ABSTRACT

An improved apparatus and method for creating high quality virtual reality panoramas is disclosed that yields dramatic improvements during the authoring and projecting cycles, with speeds up to several orders of magnitude faster than prior systems. In a preferred embodiment, a series of rectilinear images taken from a plurality of rows are pairwise registered with one another, and locally optimized using a pairwise objective function (local error function) that minimizes certain parameters in a projective transformation, using an improved iterative procedure. The local error function values for the pairwise registrations are then saved and used to construct a quadratic surface to approximate a global optimization function (global error function). The chain rule is used to avoid the direct evaluation of the global objective function, saving computation. In one embodiment concerning the blending aspect of the present invention, an improved procedure is described that relies on Laplacian and Gaussian pyramids, using a blend mask whose boundaries are determined by the grassfire transform. An improved iterative procedure is disclosed for the blending that also determines at what level of the pyramid to perform blending, and results in low frequency image components being blended over a wider region and high frequency components being blended over a narrower region. Human interaction and input is also provided to allow manual projective registration, initial calibration and feedback in the selection of photos and convergence of the system.

BACKGROUND OF THE INVENTION

[0001] 1. Field of Invention

[0002] The present invention relates generally to an improved system forcreating a full 360-degree virtual reality panorama from rectilinearimages.

[0003] A panorama is a compact representation of the environment viewedfrom a 3D position. While an ordinary image can capture only a smallportion of the environment, a panorama can capture it all, or anyportion of it, depending on the geometry in which the panoramas arerepresented. Recently there has been an explosive popularity ofpanoramas on the world wide web and in multimedia as an effective toolto present a photo-realistic virtual reality. However, creatinghigh-quality panoramas, especially those that completely enclose space,has been difficult.

[0004] 2. Description of Related Art

[0005] Various systems have been proposed for simulating a virtualreality environment using photographic quality images. Many virtualreality environments use 3D models or mathematical equations to create asimulated world. The user explores this simulation in real time. Though3D modeling via equations has certain advantages, such as a depiction ofa scene from any arbitrary vantage point, creating images from equationsgenerated by a computer is seriously limited by the speed of thecomputer. To avoid this problem, technology such as QuickTime™ VR fromApple Corporation uses images that have already been produced, eitherphotographically or generated by a 3D modeling program, and stored insecondary memory. Software only has to read the image files from a diskand display the scene as needed, rather than calculating the scene frommathematical models. However, a limitation of the QuickTime™ VR programis that it requires that the view direction for all photos reside in asingle plane, such as that obtained by rotating a camera on a tripod. Italso requires that the vertical field of view (or equivalently, thefocal length) be known, and that there be roughly equal angularincrements between one photo and the next.

[0006] Further, a panoramic movie or image can be created usingspecialized hardware, such as with a panoramic camera or a fisheye lenscamera. However, such hardware is inconvenient for the average novicephotographer. In the alternative, software can be used to simulate apanorama. This obviates the need for specialized hardware.

[0007] Though various software programs have been proposed to simulatepanoramas without the use of special hardware, these programs havecertain serious drawbacks that have not been successfully overcome todate. These include, but are not limited to, unrealistic representationsof images, lack of proper registration and calibration of images, lackof proper blending of images, and slow speed in registering, calibratingand blending images to create a panorama.

SUMMARY OF THE INVENTION

[0008] Accordingly, one aspect of the present invention is to provide animproved system and method for overcoming the drawbacks of priortechniques discussed above.

[0009] Another aspect of the present invention is to provide for theregistration, calibration and global optimization of images, preferablycaptured from a substantially single nodal position. The solution tocreating a full 360-degree panorama quickly and seamlessly is dividedinto three steps. The first step registers all overlapping imagesprojectively. A combination of a gradient-based optimization method anda correlation-based linear search has proved to be robust in cases ofdrastic exposure differences and small amount of parallax. The secondstep takes the projective matrices and their associated Hessian matricesas inputs, and calibrates the internal and external parameters of everyimage through a global optimization. The objective is to minimize theoverall image discrepancies in all overlap regions while convertingprojective matrices into camera parameters such as focal length, aspectratio, image center, 3D orientation and the like. Improved techniquesfor global optimization are disclosed that give order of magnitudeimprovements over prior systems of optimization. The third stepre-projects all images onto a panorama by a method employingLaplacian-pyramid based blending using a Gaussian blend mask generatedby the grassfire transform. The purpose of the blending is to provide asmooth transition between images and eliminate small residues ofmisalignments resulting from parallax or imperfect pairwiseregistrations. The invention further provides for human interaction,where necessary, for initialization, feedback and manual options.

[0010] Further, the present invention, unlike some of the prior art,allows for multiple views, from multiple planes and rows of images, andallows for the arbitrary orientation of photographic images to beconstructed into a panorama, without specialized hardware such as atripod or fisheye lens. In addition, the present system and method canbe several orders of magnitude faster than the prior art.

[0011] The numerous aspects of the invention described herein result ina system for registration, calibration and blending that creates highquality panoramas from rectilinear images that is up to several ordersof magnitude faster than prior systems. In one calculation, the presentinvention is up to 100,000 times faster than prior techniques. As aconsequence, the present invention could be used to construct panoramasmuch quicker than previous methods. These panoramas can be used inapplications where real-time image rendering is important, such as inreal-time 3D virtual reality, the construction of background images,computer animation, -multimedia, and the like.

[0012] The above described and many other features and attendantadvantages of the present invention will become apparent from aconsideration of the following detailed description when considered inconjunction with the accompanying drawings.

BRIEF DESCRIPTION OF THE DRAWINGS

[0013] Detailed description of preferred embodiments of the inventionwill be made with reference to the accompanying drawings.

[0014] FIGS. 1(a) and 1(b) is an artist's rendition of a compositephotograph before and after the application of the present invention.

[0015]FIG. 2 is a generalized flowchart of the various function modulescomprising one embodiment of the present invention.

[0016]FIG. 3 is a generalized flowchart of the method of operation ofthe Pairwise Registration function module of one embodiment of thepresent invention.

[0017]FIG. 4 is a generalized flowchart for the method of operation ofthe Calibration and Global Optimization function module of oneembodiment of the present invention.

[0018]FIG. 5 is a generalized flowchart for the method of operation ofthe blending function module of one embodiment of the present invention.

[0019]FIG. 6 is a screen shot of a user interface dialog window for theuser interface of one embodiment of the present invention.

[0020]FIG. 7 is a conceptual illustration on the problem of finding theproper Laplacian pyramid level using the minor axis of an inertialellipse.

[0021]FIG. 8 is a graphical illustration of the transition lengths fordifferent frequency image components (low, middle and high) used in oneembodiment of the blending function module of the present invention.

[0022]FIG. 9 conceptually illustrates the weighted average method forblending.

[0023] FIGS. 10 (a) and (b) illustrate the blend mask used for blendingtwo images during the blending phase of one embodiment of the presentinvention.

[0024] FIGS. 11(a) and (b) illustrate a particular problem overcomeduring blending in one embodiment of the present invention.

[0025]FIG. 12 illustrates a particular virtual reality orientation ofimages for the user interface for one embodiment of the presentinvention.

DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS

[0026] Disclosed herein is a detailed description of the best presentlyknown mode of carrying out the invention. This description is not to betaken in a limiting sense, but is made merely for the purpose ofillustrating the general principles of the invention. The section titlesand overall organization of the present detailed description are for thepurpose of convenience only and are not intended to limit the presentinvention.

[0027] Turning to FIG. 1, there is shown a simulation of differentoverlapping rectilinear images, or 2D photographs, framed by dashedlines during the authoring portion of the present invention, asindicated by dashed lines 110. Where the images overlap there ispotential for misalignment when constructing a 3D panorama, as indicatedby blurry lines 112, for a variety of reasons, including the arbitraryposition of the camera, errors in internal and external cameraparameters, and, distortions that occur when warping a 2D image toconstruct a 3D image space. The present invention is designed tocalibrate and align all such 2D rectilinear images with respect to oneanother and globally, blend the images where they overlap, and constructa reconstructed and relatively error free 3D panorama image, shownconceptually in 2D form as FIG. 1(b), for any arbitrary geometry.

[0028]FIG. 2 discloses a generalized flowchart for overall operation ofthe invention. The invention is part of a system 210 comprising acomputer, having all necessary hardware, such as processing unit 212,e.g., a G3 microprocessor chip; I/O 214, such as a keyboard 216, videomonitor 218, and mouse 217; memory 220, which may be any sort of memorybuffer, preferably primary memory, e.g. RAM, that can be cached tosecondary memory, such as a hard drive. The system 210 is controlled bya program residing in system memory 220, which also stores output dataand other data. The program is preferably written in the C or C++language, using classes, structures, functions, calls, translationunits, headers, subroutines, modules, and other features of structuredprogramming, where appropriate, of both data and source code, suitablycompiled and in executable form, in accordance with the teachings of thepresent disclosure to practice the invention. The present invention isalso suitable for implementation with an interpreted language such asJava.

[0029] In constructing a panorama from rectilinear images, the systemfinds solutions to three sub-problems: (1) the projective registrationsof overlapping images (shown as the “local pairwise registration” box222 in FIG. 2), (2) calibration and global optimization of these images,a self-calibration in which 2D image planes are positioned as 3D planesin space (shown as the “calibration and global optimization” box 224 inFIG. 2), and (3) the composing or blending problem in which images areready to be reprojected to a 3D environment map with pixels in overlapregions being composed from multiple images, to smooth any transitionaldiscontinuities (shown as the “blending” box 226 in FIG. 2). Finally,there is the projection or construction of the assembled panorama onto a3D geometry surface, such as a cylinder, cube or sphere (defined as the“projection” box 228 in FIG. 2).

[0030] The solutions to these sub-problems are performed by softwarefunction modules 222, 224, 226, 228 residing in memory 220 and operatingthe processor 212. The modules are designated, as explained furtherherein, the pairwise registration function module 222, the calibrationand global optimization function module 224, the blending functionmodule 226, and the projection function module 228. A user interfacemodule 230, also residing in memory 220, may interact with the othermodules to pass data to and from the modules, and accept input from ahuman user of the system. The modules may receive data from memory,manipulate that data as described herein, and output the data to othermodules. The three modules 222, 224 and 226 may perform feedback to passdata back to previous modules, as indicated by arrows 233, and asdescribed below. Although in a preferred embodiment the modules areprogrammed as separate software routines or classes, the modules may becombined into one module performing all the designated tasks performedby separate modules.

[0031] As a final step, the fourth module, the projection functionmodule 228, constructs a panoramic scene by projecting the blended imageonto any designated geometry view surface, typically a cubic,polyhedral, cylindrical or spherical surface. The projection module maybe controlled through the user interface 230 as well, to allow a user toselect what geometry will be projected onto and to control and modifyother factors, including the use of photo re-touching software such asPhotoShop™ for modifying the final panorama.

[0032] Generally, the local registration, self-calibration and globaloptimization, and blending involve a multi-step procedure.

[0033] First, regarding the initial local registration, and referringgenerally to the generalized flowchart of FIG. 3, the system 210 readsin each overlapping rectilinear image into main memory 220, as indicatedby step 312. The images are assumed to roughly share a common nodalpoint (i.e., that point in the three-space where all rays of lightconverge through a lens) with other overlapping rectilinear images. Theobject of the program during local registration is to register thelocally overlapping images, by comparing common overlapping areasbetween overlapping images at certain predetermined resolution levels-on a Gaussian pyramid representing the overlapped images. Differentcombinations of overlapping areas are tried to achieve the optimaloverlap between images (or, equivalently, the smallest error in theerror function or pairwise objective function described herein) usingthe steps described herein, which generally minimizes the averagesquared pixel intensity (e.g., brightness and contrast) difference withrespect to certain transformation parameters. Initial values forparameters used in optimizing the pairwise objective function areassumed by the computer, as indicated in step 314. The initial valuesmay optionally be input by a user, e.g., with a user interface 230 as inFIG. 2, and in response to a user dialog window such as of the kindshown in FIG. 6. Besides the global orientation (pan, tilt and roll) theother parameters that are most likely to give instability in convergenceof the error function are bad initial estimates of the brightness andcontrast, as well as of the geometric image center of projection of theoverlapping rectilinear images. Certain parameters most likely to createinstability in the convergence of the local error functions can becontrolled (e.g., progressively dampened at different levels of theGaussian pyramid) to ensure convergence, as indicated by step 318. Theoverlapping images are then perturbed and the local error function withrespect to these and other variables is calculated until a minimal localerror function is found, as indicated by step 320. The minimal localerror function is then stored for a particular level of the Gaussianpyramid, as indicated in step 322, for each pairwise registration, andis saved and later used to compute a global error function for all theoverlapping images. The local pairwise registration module 222 iteratesuntil the entire Gaussian pyramid is traversed, starting from thecoarsest level of the pyramid (sometimes called the bottom, where thepyramid can be standing on its inverted top) and working to the finestlevel resolution, as indicated in decision box 324. It should be notedthat at any stage throughout the registration, and throughout theinvention in general, the system may check for a user interruption,through the user interface, that would require immediate attention fromthe processor, such as to allow the user to interactively adjust theparameters to avoid divergence or convergence to an undesired localminimum.

[0034] As indicated in box 316 of FIG. 3, it must be determined at whatlevel in the Gaussian pyramid to start the local pairwise registration.One way to find the lowest level is to select the resolution level atwhich it is found that the images share at least some arbitrary numberof overlapping pixels, e.g., preferably about 30 pixels from each side,e.g., preferably no less than 30 pixels across the overlapping area. Ifgreater than 60 pixels of overlap is found in these areas, the size(resolution) of the overlap region is decreased by half (going deeperinto the pyramid) and the procedure of the present invention isreiterated again. If, on the other hand, the overlap is less than 30pixels, then the size (resolution) of the overlap region is increased bydoubling. By utilizing multi-resolution registration of overlappingimages by way of the Gaussian pyramid, convergence to the desiredoptimum is accelerated, and false local minima are avoided.

[0035] On occasion, it may be visually apparent to a user that duringregistration the images are not converging optimally. In this case userinput may manually abort the pairwise registration procedure, and theuser may manually help align the images closer before resuming automaticregistration, as before. This manual intervention is true for allaspects of the invention. Nevertheless, the present invention issurprisingly robust, and manual intervention is not a prerequisite forthe invention to work.

[0036] Non-optimal convergence or divergence has sometimes been found tobe the case whenever images for a spherical projection are used,especially those in the “pole” regions of the sphere (though in generalthe invention can adjust quite nicely for images that wrap around thepoles). Divergence sometimes results when the initial default parameterschosen are wildly off or not suitable for convergence. During suchinstability, the images will appear to a user to “run away” from eachother. In this case, and throughout the invention, provision may beprovided in the user interface 230 of the embodiment of FIG. 2 of thepresent invention for manual intervention, such as to abort the program,for the manual selection and relative positioning of the images to bepairwise registered, and for the selection and relative positioning ofoverlapping images for blending.

[0037] The iterative method of moving down a pyramid when an overlapregion is greater than, say, 30 pixels, is an attempt to preventinstability in the error function due to problematic parameters, such asinitial value errors in the image center of projection of the imagesbeing registered, and errors in setting initial brightness and contrastvalues. Techniques of damping and annealing of problematic parameters(with damping progressively diminished and finally set to zero as onemoves up the pyramid to finer levels) can be used to stabilize the localerror function for these problematic terms, as explained further herein.

[0038] One improvement over prior techniques has been to save the localerror function values and use them to compute and optimize the globalerror function needed for optimization. This improvement also avoidshaving to evaluate the entire global error function (global objectivefunction) from scratch. The pairwise objective functions (local errorfunctions) are approximated by a quadratic Taylor series, and, togetherwith the chain rule, the global objective function (global errorfunction) is minimized. Calculation of the global error function isgreatly speeded up by this procedure.

[0039] Further regarding rectilinear images taken in a non-arbitrarymanner (such as from a tripod that is rotated, or a photographer whomanually “pans” a field of view), the number of pyramid levels andoptimal direction for the blending of images overlapping in a region canbe computed by the present invention by computing the minimum eigenvalueof the 2×2 inertial tensor of the overlapping region between two images.It has been found in practice that for an arbitrary polygon representingan overlapping region, the optimal direction for blending, as well asthe width of the blending region (which determines the level in thepyramid at which to start the method of registration and optimization)is found along the minor axis of the inertial ellipse found from solvingfor the inertial tensor of the overlapping images. A similar method offinding the proper pyramid level is by solving for the smallesteigenvalue of an inertial tensor of the overlap region between images.Conceptually, such an ellipse is shown in FIG. 7. The blending region inan arbitrarily shaped polygon region 710, which represents the area ofoverlap between overlapping images, lies along the width and directionof the minor axis 712 of the ellipse 714, which is calculated from theinertial tensor of the overlapping images forming the polygon.

[0040] Thus, the results from computing the inertial tensor are used todetermine the pyramid level, blending width, and blending direction. Thesmallest inertial eigenvalue is used to determine the number of pyramidlevels. One could also use the eigenvalue vector (eigenvector) todetermine the direction, or, preferably, use a blending mask, asexplained herein, that yields a grayscale ramp, which defines directionin a direction field from taking the grayscale ramp gradient.

[0041] Next, after pairwise local registration, global optimization isused to remove any inconsistencies. The parameters found at the localregistration level, generally from six to ten parameters per overlappingimage pair, are optimized globally. Various constraints ensureoptimization whenever there are suitable overlapping image pairs, as thenumber of independent parameters is usually less than the constraints.Regarding optimization in general, to ensure the best chance forconvergence of a solution a combination of simulated annealing andprogressive damping is used, as described herein.

[0042] Global optimization is necessary because noise in the images willyield inconsistencies in cyclically overlapping sets of images (e.g.,that A is pairwise registered with B, B with C, and C with A, does notnecessarily mean A and C are property globally registered with C and B):During the global optimization phase, the discrepancies are distributedamong all image pairs in such a way as to minimize the increase in totalerror. The Hessian matrix, computed in the local registration phase, asdescribed further herein, provides a description of the shape of theerror function landscape in the neighborhood of each pairwise optimum.Parametric perturbations occur along the valleys of this landscape (asindicated by the Hessian), where the increase in error is minimal. Oneway of looking at this solution is to say that knowledge gathered fromthe pairwise registration optimization is saved and used for globaloptimization, avoiding additional computation.

[0043] Turning attention to FIG. 4, there is shown a generalizedflowchart for the calibration and global optimization module of thepresent invention. Data from the local pairwise registration module 222is provided to the calibration and global optimization module 224, asindicated by box step 402. Such data can include the pairwise objectivefunctions (local error function values) for each pairwise registrationfound previously. The global objective function is calculated from suchdata, as indicated by boxes 404, 406, and as described more fully below.The alignment of images globally is checked to ensure globalregistration, as in decision boxes 410 and 412, which may accept manualinput from a user, via a user interface module. If there is alignment,the system proceeds from the calibration and global optimizationfunction module to the blending of images. Otherwise, as illustrated byboxes 418 and 420, the pairwise registration module parameters (such asshown by box 314 in FIG. 3) may be re-initialized and the pairwiseregistration module reexecuted to recompute the pairwise registration ofimages, using better, updated camera parameters as determined from theglobal registration procedure as shown in FIG. 4.

[0044] During the blending step of the procedure utilized by the systemof the present invention, the image overlap regions are “blended” or theimages are “stitched” together, so that the high frequencies (e.g.,sharp lines of contrast, such as edge boundaries, analogous to the highfrequency signals associated with a square wave) are blended over anarrow blend region range, and the low frequencies (e.g., illuminationvariations, analogous to DC baseband signals) are blended over a widerange. In this way the images are seamlessly integrated to form apanorama in an aesthetic manner.

[0045] In one preferred embodiment of the invention, as illustrated bythe steps in FIG. 5, blending is performed by determining the coarsestlevel in a Laplacian pyramid for which to begin blending of images I, J(box 504), and constructing a Laplacian pyramid at this level (box 506).One of the two images may be comprised of previously blended images. Ablend mask boundary is generated (box 508), defining the boundary overwhich blending is to occur, preferably using the “grassfire” transformmethod. Next, a blend mask is generated by a Gaussian pyramid method(box 508). In one preferred embodiment, the overlapping images I and Jthat are to be blended are put into a Laplacian pyramid and multipliedby the blending mask, or its compliment, respectively (box 510). Theresulting product is repeatedly added together to each level of theLaplacian pyramid (box 512), moving up to the finest level resolution ofthe pyramid in a sequential fashion (decision box 514), until a blendedimage is achieved for the two images. A similar procedure is performedfor all other images that overlap, as indicated by decision box 516.This preferred technique has the advantage over prior techniques in thatlow frequency component images are blended over a wider area region,giving a smoothing effect, as desired for low frequency components,while high frequency components (such as sharp edges) are blended over asmaller blend region, giving a “sharpening” effect for these highfrequency components, exactly as desired. This is illustratedconceptually by FIG. 8.

[0046] Authoring a panorama from 2D images can be thought of as dividedinto two different phases: (1) orientation of originally 2D images into3D space, and (2) the projection of a panorama onto a particular 3Dgeometry, that can later be used to project views of the panorama onto a2D viewing plane. A series of preferably overlapping photographs areanalyzed to determine what orientation the photographs were taken inorder to establish a common ground for subsequent operations, includingthe construction of a panorama. The panorama is not ordinarily meant tobe viewed by a user, only the subsequent projection of the panorama ontoa viewing plane is viewed by the user. The panorama is constructed on aparticular geometry that will best facilitate the subsequent step(sometimes termed rendering) of the projection of the panorama from theparticular geometry onto a chosen viewing plane for viewing by a user.Typical geometries in the system of the present invention on whichpanaoramas are formed include: cubic, polyhedral, cylindrical andspherical geometries. However, any type of geometry may be used, such astwo frusto-conical cones joined at their base with the apexes pointingaway from one another; any quadric surface, and any and all of thegeometries that employ the following projections: equidistant,equiangular, ellipsoid, Mercator (and all derivatives thereof, e.g.,transverse Mercator, Oblique Mercator, and the like), cylindricalequal-area, Miller cylindrical, equidistant cylindrical, Cassini (e.g.,both for spherical and ellipsoid projections, and the like), all conicmap projections, e.g., Albers equal-area, Lambert conformal conic,equidistant conic, bipolar oblique conic conformal, polyconic, Bonne,all azimuthal and related projections, e.g., orthographic, sterographic,gnomonic, general perspective, Lambert azimuthal equal-area, azimuthalequidistant, modified-stereographic conformal, all space mapprojections, including space oblique Mercator and satellite-trackingprojections, all pseudocylindrical and other miscellaneous projections,including Van der Grinten, sinusoidal, Mollweide and Eckert IV and VIprojections. The foregoing list is meant to be illustrative and notexhaustive of the geometries and projections possible during theconstruction and employment of panoramas using the system of the presentinvention.

[0047] Further, the present invention can be employed in future systemsthat are fast enough to eliminate the need for a projection functionmodule 228, and proceed directly from pairwise registration, calibrationand global optimization to the viewing and blending of the panorama on achosen viewing plane, without loss of generality. Presently, however,there is not sufficient computing power in most desktop computers forthis to be feasible for real time applications.

[0048] Furthermore, the system of the present invention has means for auser interface for all phases of the invention. A user may select, amongother things, which images are to be registered, and at what arbitraryimage plane. The user interface, suitable for display on a computermonitor and with input from a keyboard, mouse pointer, or other I/Odevice, has fields for any and all internal and external parameters ofthe projection matrix of the images, including aspect ratio, number ofrows of images, the tilt between rows, the angle between photos within arow, the roll of each image taken (e.g., landscape mode), as well asfields for how many horizontal rows of images are to be registered(typically two or more), image center position, focal length of camera,camera orientation with respect to a common reference frame, such ascamera pan, tilt, roll and skew, and the brightness and contrast ofimages. The user interface may have the ability to adjust theaforementioned parameters for each image individually, or may have theability to adjust parameters for images captured with a particularmethodology, such as equal angular increments in latitude and longitude.

[0049] Thus, turning attention now to FIG. 6, there is shown a screenshot of a user interface for the present system, suitable forinitialization parameters for an ensemble of images, captured with anequal angular increment methodology that can be facilitated by use of atripod. The user interface is particularly useful for authoringpanoramas when a user wishes to adjust the automatic default orientationgenerated by the computer for the present invention. A more complex userinterface may be provided to accommodate the initial orientation andplacement of a free-form set of images, i.e. a set of images capturedwithout any particular methodology, such as that captured using ahand-held camera. A plurality of parameters of the kind described hereinmay be manually entered into the dialog box fields by the user (ifknown) to aid in the pairwise registration, calibration and globaloptimization and blending of images. Further, the particular imagesselected for pairwise registration and calibration and globaloptimization may be aborted during non-convergence run-away conditions.Some of the parameters that may be explicitly specified by a userinclude (referring to FIG. 6) the number of rows 602 (assuming a panningof photos are taken, with a number of overlapping rows of photos takenabout a 360 degree arc), initial pan 604, initial tilt 606, panincrement 608 and initial roll 610 for the first row of photos (rows arepreferably used, but vertical columns of photos are also contemplated);and, for inter-row parameters, the pan increment 612, roll increment 614and tilt increment 616. Zoom lens distortion factors and othermiscellaneous factors may be entered in dialog fields such as field 620.Other parameters may be specified for overriding the computer defaultsand for better guaranteeing convergence, such as camera focal length orf-stop, pixel aspect ratio, and the like. In addition, the userinterface may allow for the selection, arrangement and relativepositioning of photos to be composed into a panorama, preferably in arow by row layout, with preferably at least one or more rows of photosto be made into the panorama.

[0050] Further, while one preferred embodiment of the present system isdesigned for overlapping image pairs to share a common nodal position,in general the system is forgiving of nodal point changes, provided thechanges are not excessive, e.g., a 1% shift in nodal point should stillallow satisfactory panoramas to be formed. The one-percent shift is notoverly cumbersome when taking long range photos, e.g., from the heightof the 300 m Eiffel Tower a 1% shift will allow up to 3 meters (nearly10 ft) of shifting of camera nodal position when taking photos at thebase, which can accommodate an amateur photographer taking pictureswithout a tripod.

[0051] Turning attention again to the three modules labeled PairwiseRegistration, Calibration and Global Optimization and Blending, asillustrated in FIG. 2, the authoring aspect of the invention will befurther described.

[0052] I. Pairwise Registration

[0053] To find a solution to the first sub-problem posed in constructinga panorama from rectilinear images, finding the projective registrationsof overlapping images, one must pairwise register the two images.Pairwise registration can be thought of as synonymous to finding anestimate of the projective transformation relating two given overlappingrectilinear images. The projective transformation may be represented bya particular parametrized projective matrix, that is parametrized by atypical canonical number (usually 8 or 9) projective parameters, e.g.,3D rotation parameters (pan, tilt roll), center of projection of images,ratio of focal lengths, and the like. The projective matrix can bedefined as a particular case of a three-dimensional affinetransformation, a transformation that effects rotation, scaling, shearand translation, with the restriction that camera motions arerotational. As is known per se, transformations are important tools ingenerating three-dimensional scenes, in moving objects around in anenvironment, and in constructing a two-dimensional view of theenvironment.

[0054] The mathematics described herein are but one representation forthe method carried out by the apparatus of the present system for one ormore preferred embodiments of the invention. In mathematics the samephenomena can be represented in different symbolic notation—which oftenappear to the untrained eye as being radically different from oneanother—without changing the nature of the phenomena described. Forexample, as explained above and represented further below, theprojective transformation of a pair of images can be more particularlycharacterized as a projective transformation that is particularized as aprojective matrix parametrized by a certain projective matrix parameters(typically having 8 or 9 projective parameters, as explained herein,such as pan, tilt, roll, center of projection of the images, ratio offocal lengths, and the like). However, this particular representationdoes not preclude the projective transformation from being reduced topractice using the teachings of the present invention by alternateequivalent methods or other representations, other than as representedby a particular parametric matrix representation, without loss ofgenerality from the way the invention is described herein. Further, andconcomitantly, the transformations involved with the present inventionmay be described in alternate notation, using for example Euler anglesor quaternions, without detracting from the spirit and scope of theinvention. It is to be understood from the teachings of the presentdisclosure that the description of a projective matrix also includesthese other representations. By the same token, programming constructssuch as data structures and classes are typically realized in binarycode, rather than abstract mathematical notations, and, as such,constitute the machine readable representations of the constructs. Therepresentation of such constructs in this form do not result in any lossof generality of the representation of the underlying invention asdescribed herein.

[0055] Regarding local pairwise registration in general, if onerestricts camera motions to be rotational only, the 2D warping betweenimages i, j, is strictly projective in absence of lens distortions, andgiven by, [Eq. (1)] $\begin{matrix}{\begin{bmatrix}x_{i} \\y_{i} \\z_{i}\end{bmatrix} = {\begin{bmatrix}m_{0} & m_{1} & m_{2} \\m_{3} & m_{4} & m_{5} \\m_{6} & m_{7} & m_{8}\end{bmatrix}\begin{bmatrix}x_{j} \\y_{j} \\z_{j}\end{bmatrix}}} & (1)\end{matrix}$

[0056] where:[_(X) _(i) _(y) _(i) _(Z) _(i) ]^(T) are the homogeneouscoordinates of pixel locations (with the convention that column vectorsrepresent three-dimensional points). In the following description, thevector X_(i) represents the homogeneous coordinates, and the matrixM_(ij) represents the matrix that transforms coordinates from image j toimage i. Due to the scale ambiguity in the projective matrix, the lastparameter m₈ in the projective matrices is set to equal 1.

[0057] The objective of local pairwise registration is to estimate theprojective matrix given two overlapping images. The projective matrix isinitialized by the camera internal and external parameters, e.g.,[Eq.(2)]

M _(ij) =T ⁻¹(p _(i) q _(i))T(p _(j) ,q _(j))  (2)

[0058] [Eq. (3)]

[0059] where $\begin{matrix}{{{where}\quad {T( {p_{i},q_{i}} )}} = {{R( q_{i} )}\begin{bmatrix}1 & 0 & {- C_{x}^{i}} \\0 & a_{i} & {- C_{y}^{i}} \\0 & 0 & f_{i}\end{bmatrix}}} & (3)\end{matrix}$

[0060] where [C_(x) ^(i), C_(y) ^(i)],a_(i),f_(i) are the image centerposition, the aspect ratio and the focal length, respectively;

p _(i) =[a _(i) , f _(i) , C _(x) ^(i) , C _(y) ^(i)]^(T)

[0061] is the internal parameters, vector;

[0062] q_(i)=represents the camera orientation. with respect to a commonreference frame; and

[0063] R( )=represents the 3×3 rotation matrix computed from theorientation parameters q_(i).

[0064] Camera internal and external parameters are initialized eitherautomatically by the computer assuming default values, or manually withuser input.

[0065] There are ten parameters in the projective registration: eightindependent parameters in the projective matrix and two parameters tocompensate for brightness and contrast difference between the twoimages. The gradient-based optimization minimizes the followingobjective, as suggested in box 320 in FIG. 3, by instructing theprocessor 212 to perturb the overlapping images stored in memory withvarious combinations of overlapping pixels until the below localregistration error function has the smallest value [Eq. 4]:$\begin{matrix}{e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}} & (4)\end{matrix}$

[0066] where s_(ij) and b_(ij), the exposure parameters, represent theexposure difference, I_(i)( ) and I_(j)( ) are pixel intensity valuesfrom the two images, and A_(ij) is the overlap area (which helpsnormalize the error function e_(ij)). The optimizations are performed onprogressively finer levels of Gaussian pyramids. In practice, however,it has been found that the direct application of gradient-basedoptimization frequently failed due to exposure differences, largetranslations, or both. Therefore, preferably a combination ofcorrelation-based linear search and a progressive damping (e.g.,simulated annealing) of exposure parameters is used to alleviate theproblem, as suggested by box 318 of FIG. 3. On the coarsest pyramidlevel (e.g., ‘lowest’ resolution level of the image, the method of theinvention first performs a linear search over the transitionalparameters using normalized correlations, an idea similar to aprogressive complexity search known per se in the art, e.g., see H.Sawhney and R. Kumar, “True multi-image alignment and its application tomosaicing and lens distortion correction”, Proc. of CVPR, pp. 450-56(1997). Since the image size on the coarsest pyramid level is small, thecorrelations are done efficiently. Once the maximal correlations arefound, the exposure parameters s_(ij) and b_(ij) are estimated through alinear regression. When the gradient-based optimization is performed onsubsequent finer pyramid levels, the damping coefficients on exposureparameters are reduced exponentially, and finally set to zero at thefinest pyramid level.

[0067] To determine the number of pyramid levels given an arbitraryoverlap of two images, one can compute the eigenvalues of the 2×2inertial tensor of the overlap polygon region. Determining eigenvaluesand inertial tensors are known in the art per se. The square root “I” ofthe smaller eigenvalue is used to estimate the number of pyramid levelsaccording to the formula: $\begin{matrix}{\log_{2}( \frac{l}{l_{\min}} )} & (5)\end{matrix}$

[0068] where I_(min) is the minimal size of the finest level resolutionpyramid level. In a preferred embodiment, I_(min) is set to 10 pixels.

[0069] II. Calibration and Global Optimization

[0070] The second major step in authoring panoramas is to extract camerainternal and external parameters from those projective matricesestimated in step I. above. In general, it is impossible to achieve adirect solution by inverting Eq. (2) above directly to obtain the cameraparameters, since there are eleven camera parameters while a projectivematrix provides only eight constraints. However, because one imageusually overlaps with multiple images, one can take advantage ofredundancy in the system to obtain a consistent set of camera parametersthat approximates all projective matrices in the same time. A globaloptimization routine module is used to achieve this goal.

[0071] Since the projective matrix is a function of camera parameters asin Eq. (2), to extract all camera internal and external parameters thefollowing objective functions are minimized by having the calibrationand global optimization module 224 of FIG. 2 instruct the processor toperturb different combinations of images stored in memory to minimizethe following global error function: [Eq. (6)]

E=Σ _(ij) A _(ij) e _(ij)(M _(ij)(p _(i) ,q _(i) ,p _(j) ,q _(j)))  (6)

[0072] where e_(ij) is the pairwise objective function in Eq. (4).

[0073] However, it has been found that it is computationallyprohibitively expensive to evaluate the objective functions according toEq. (6). By noting that the pairwise objective function e_(ij) hasalready been optimized individually, we can approximate it by aquadratic surface, which can be viewed as terms from a matrix TaylorSeries expansion. Thus, in a preferred embodiment, the followingapproximation is used in the calibration and global optimization routinemodule: [Eq. (7)]

e _(ij)(M _(ij))≈e _(ij) ⁰+(M _(ij) −M _(ij) ⁰)^(T) C _(ij)(M _(ij) −M_(ij) ⁰)  (7)

[0074] where e_(ij) ⁰ is a constant representing the minimal valueachieved in the pairwise registration; M_(ij) ⁰ is the 8×1 vectorrepresenting the optimal projective matrix, and

[0075] C_(ij) is the 8×8 Hessian matrix obtained when optimizingobjective function e_(ij), as in the methods of Levenberg-Marquardt orBroyden-Fietcher-Goldfarb-Shanno, known in the art per se. Thus theprojective matrix can be treated as an 8×1 vector instead of a 3×3matrix. This step is represented in FIG. 4 as box 404 Other methods maybe employed in the present invention using the teachings herein withoutdeparting from the scope of the invention.

[0076] Next, once the pairwise objective functions are approximated byquadratic surfaces with Eq. (7), the global objective function in Eq.(6) is used as a weighted sum of all those quadratic surfaces thuscomputed in the global optimization routine module. The global objectivefunction has a gradient with respect to the camera internal and externalparameters that can be easily established thorough the chain rule: [Eq.(8)] $\begin{matrix}{\frac{\partial E}{\partial( {p_{i},q_{i}} )} = {\sum\limits_{j}{\frac{\partial e_{ij}}{\partial M_{ij}}\frac{\partial M_{ij}}{\partial( {p_{i},q_{i}} )}}}} & (8)\end{matrix}$

[0077] from using Eq. (2) and Eq. (7) above, as shown in FIG. 4 as box406. It has been found that since no direct evaluation on images isinvolved, the computation required using Eq. (7)-(8) as described abovein minimizing the global objective function for all the overlappingimages is nearly trivial. The computational savings using the techniquesdescribed herein result in several orders of magnitude in savings oftime and speed over prior techniques of optimization.

[0078] In the most general case, the camera parameters for each imageare four internal parameters, pi, and three orientation parameters qi.Every pairwise registration provides eight constraints on those cameraparameters. When there are plenty of overlapping images pairs, theoptimization is overconstrained, in that the number of independentparameters is generally less than that of the number of constraints.However, in practice, even when the optimization appears to beoverconstrained, oftentimes camera parameters are so weakly constrainedthat they can easily diverge the whole optimization.

[0079] In order to solve the aforementioned problem, so that theoptimization described herein behaves well in underconstrained or weaklyconstrained situations, it is preferred that simulated annealing be usedto dampen the camera internal parameters. Conceptually this step isshown in the generalized flowchart of FIG. 4 as box 408. As theoptimization progresses, the damping parameters are gradually reduced.The exposure parameters s_(ij) and b_(ij) are estimated through a linearregression. The proposed solution has been found to work remarkably wellin practice, and is a significant improvement over prior techniques.

[0080] The pairwise registration and global optimization described abovecan be iterated if the alignments are still not satisfactory, by eitherby a user acting through the user interface, or by automated means, asindicated by arrows 233 in FIG. 2 and boxes 418 and 420 in FIG. 4. Thus,in the iteration, the pairwise registration module 222 will use theimproved camera parameters generated by the calibration and globaloptimization module 224 to re-initialize the projective registrations,and re-compute the optimal projective matrices and their Hessians. Thenthe improved projective matrix parameters will, in turn, be used togenerate improved estimations of camera parameters in the globaloptimization module 224. Likewise, blending may be further iteratedafter registering, calibrating and optimizing the images after aninitial blending. III. Blending

[0081] Notwithstanding the improved method and system of pairwiseregistration, camera calibration and global optimization of imagesdescribed above, it has been found that for high quality panoramasblending of overlapping image boundaries is often required prior toviewing the panorama. Such blending is required when the pairwiseregistration and global optimization by the system generates panoramashaving imperfectly aligned images that give “shadow” or “ghosting”effects, if the images are averaged in overlap regions. As human eyesare very sensitive to such imperfections, in one preferred embodimenthuman input may be used to facilitate proper blending during theauthoring of a panorama, using the system and method disclosed herein.Other methods for blending may also be employed in the presentinvention, such as the multi-resolution weighted average method, and theweighted average method. In the weighted average method, as illustratedin FIG. 9, there is a transition region 902 between images to beblended, Image 1 and Image 2. The weights of image 1 for intensity orother parameter is linearly decreased from a value of 1.0 to 0.0 in thetransition region, while the weights of Image 2 is increased from 0.0 to1.0. A pixel in the transition area is a weighted sum of two pixels fromtwo images.

[0082] By contrast, the multi-resolution weighted average method firstdecomposes two images into different frequency bands by buildingLaplacian pyramids, and performs separate weighted averages on eachpyramid level with different transition lengths for each frequency.Transition lengths are defined by the region it takes for a parameter togo from value 1 to 0 or 0 to 1. FIG. 8 shows the transition lengths fordifferent frequency bands (low, middle and high) of images, with highfrequency image components having a shorter transition length regionthan low frequency image components. The result of this multi-resolutionblending method is seamless and absent of shadow effects.Multi-resolution weighted average blending, known per se in the art, isfurther described in P. Burt and E. Adelson, “A multiresolution splinewith application to image mosaics”, ACM Transactions on Graphics,2(4):217-236 (1983). While multi-resolution blending described herein isthe preferred blending technique in a preferred embodiment of theinvention, other types of blending, including simple weighted averageblending, is within the scope of the invention.

[0083] To determine the boundary of overlap regions to performmulti-resolution blending, preferably a blend mask is needed for anarbitrarily shaped transition region. The Gaussian pyramid of the maskimage supplies the weights for every pixel at every pyramid level. FIG.10 (a) and (b) illustrate the blend mask 1000 used for the panoramiccanvas having two overlapping images, Image 1002 (Image 1) and 1004(Image 2). In order to maximize the size of the transition region forblending, the boundary curve of the mask inside the overlap regions,boundary 1001, needs to be as far way as possible from the originalimage boundaries. To locate the mask boundary a grassfire transform ispreferably used on two images individually. The resulting distance mapsrepresent how far away each pixel is from its nearest boundary. Thepixel values of the blend mask is then set to either 0 or 1 by comparingthe distance values at each pixel in the two distance maps.

[0084] Generally, the grassfire transform measures distance by notingthe time it takes for a constant velocity wave front to reach a certainpoint, knowing that distance equals velocity multiplied by time.Conceptually, the distance is measured as if distance were measured bynoting the time elapsed that a grass fire, having constant velocity andlit on the edges of a grass field, progresses to the center of thefield. In this way the exact boundaries of the field do not have to besurveyed by more precise techniques of geometry. Further detailsregarding the grassfire transform, known per se in the art, can be foundin the literature, e.g., C. Arcelli, L. P. Cordella, and S. Levialdi, “Agrassfire transformation for binary digital pictures.”, ICPR74, pp.152-54 (1974). However, though the grassfire transform is preferred inone embodiment of the present invention, other techniques may be used togenerate a blend mask, such as solving directly a Laplacian differentialequation with predetermined boundary conditions to directly generate amask with the addition of a gray scale.

[0085] Regarding the blending procedure, the blending is achieved by thefollowing method. An empty panoramic canvas, which can simply be abuffer of memory in a computer system, such as memory 220 in FIG. 2, iscopied with the first image, Image 1. Second, new images are blendedonto the panoramic canvas one by one. For each of those new images, theblend mask from the panoramic canvas and the new image is generated.Next, there are computed Laplacian pyramids of the images and Gaussianpyramids of the masks in the bounded rectangular areas of the overlapregions. Multi-resolution blending is used to blend the Laplacian andGaussian pyramids. Finally, the blended images is copied onto thepanoramic canvas, which may hold other images.

[0086] The flowchart for this procedure, generally speaking, is alongthe lines as shown conceptually in FIG. 5, and as can be modified fromthe teachings of the present invention.

[0087] (1) determine the coarsest level resolution level of a Laplacianpyramid at which two images I_(i), I_(j) (or I, J) are to be blended(step 504). The bottom-most level is computed from computing theinertial tensor of the images I, J, which gives an inertial ellipsehaving a minor axis that is used to find the number of pyramid levels inthe Laplacian pyramid.

[0088] (2) From the overlapping images I, J construct a Laplacianpyramid (step 506), as is known per se in the art, starting at thefinest level of resolution and working through the pyramid to morecoarser levels of resolution.

[0089] (3) Generate a blend mask (step 508), preferably using thegrassfire transform described herein, and construct a Gaussian pyramidfrom the blend mask. The Gaussian pyramid may be constructed by applyinga low -pass filter to the blend mask, which dilutes the sharp edges,from linear interpolation between the black and white regions of theblend mask, or from other techniques.

[0090] (4) At each level in the Laplacian pyramid of images I, J thatare to be blended, blend the overlap regions, preferably by the methodof multiplying the Laplacian values of images I, J times the weightedvalues supplied by the Gaussian pyramid based blend mask, according tothe value for I times the mask value (MV), or the value for J times thecompliment of the mask value (1−MV), (step 510).

[0091] (5) Add the results at each level of the Laplacian pyramid (step512), until the finest level resolution of the Laplacian pyramid isreached (or upon user interruption, as always) (steps 512 and 514).

[0092] (6) Perform steps (1)-(5) for all images to be blended in theblend region (step 516). Previously blended images may be blended withnew images, and blending may be iterated, with human input orautomatically.

[0093] The net result of these steps is that the lowest frequencycomponent images are blended over the entire transition region, giving adesired smoother effect, while high frequency components (such as sharpedges) are blended over a smaller region, typically 3-5 pixels wide,which “sharpens” these high frequency components, as desired.Conceptually this is shown in FIG. 8, with the blending regions 802,804, 806 for the low, middle and high frequencies being shorter,respectively, according to the functioning of the blending module asdescribed.

[0094] The net visual effect of the blending described herein, besidesproducing a pleasing smoothing of images that helps eliminate ghosting,is an “auto-iris” effect, in that luminosity appears to be automaticallyadjusted when moving from a high intensity image to a low intensityimage, similar to the adjustments made by a human eye.

[0095] In practice, it has been found that on occasion, where there areoverlaps of more than two images, a particular problem is encountered inthe blending of images, as illustrated graphically in FIG. 11(a). When afirst (N−1) images are sequentially blended onto the panorama canvas,and there is an attempt to blend Image N, most of the area covered byImage N is already blended by an “unintended” overlap between Image 1and Image N−1, with the blend mask 1110 indicated generally as the grayarea inside the dotted lines. As a result, Image N has little effect onthe panorama even though it provides much large transition areas betweenImage 1 and Image N−1, and therefore has the potential to improve thequality of the panorama. In the most general case, because of the natureof panoramas, the problem illustrated by FIG. 11(a) has the potential tooccur on occasion.

[0096] One method to solve the problem illustrated by FIG. 11(a) duringauthoring of the panorama is to allow manual input to override thecomputer default ordering for the blending of images. Images with largeroverlaps (e.g., a very dark, low frequency image that may be abackground color to a series of high frequency lighter foregroundimages) should be blended onto the panoramic canvas first manually, andcan be selected for blending by a user of the present apparatus, througha suitable user interface that lists the images to be blended first in apriority list.

[0097] Another more automated method to solve the problem illustrated byFIG. 11(a) is illustrated conceptually in FIG. 11(b). The automatedsolution uses a labeling scheme. For every pixel on the panoramiccanvas, which may be a blended pixel, the computer labels it with anumber indicating which source image contributes the most to this pixel.In FIG. 11(b), the dashed line 1120 represents the blend mask boundarywhen Image 1 and Image N−1 are blended. The pixels on the left of thedashed line 1120 have label 1, as Image 1 contributes the most to thepixel values in that area of the panoramic canvas, while the pixels onthe right side of the dashed line 1120 have label N−1, as image N−1,which may be a previously blended image, contributes the most to thepixel values in that area. When another image, Image N, needs to beblended onto the panoramic canvas, first a grassfire transform on thepanoramic canvas is performed. In addition to the actual imageboundaries, the dashed line 1120, representing the blend mask boundarybetween images Image I and Image N−1, is regarded as a virtual boundary,or “firewall”, that the grassfire cannot penetrate. The virtualboundaries are computed easily using the pixel labels and the list ofall intended overlaps. The resulting blend mask is illustrated as thegray area 1130 in FIG. 11(b), which is much larger than the blend maskarea of FIG. 11(a). Using this “firewall” technique, the blending takesadvantage of both large overlaps between Image I and Image N, and ImageN and Image N−1, which alleviates the aforementioned problem associatedwith FIG. 11(a).

[0098] IV. HUMAN INTERFACE

[0099] Human interaction is an integral part of the method and apparatusof the present invention. The system of the present invention is tryingto solve a complicated nonlinear optimization problem. No automatedprocedure can guarantee its convergence to a global minimum; in thissense the system is heuristic. On occasion, human interaction through auser interface can steer the system to produce a more satisfactorypanorama. Some areas of possible human interaction with the systeminclude manual projective registration, initial calibration, andfeedback generally.

[0100] In manual projective registration, human interaction is sometimesbeneficial in cases where the projective registration of images breaksdown due to factors such as excessive exposure difference, motion in thescene, bad initial estimates, and the like. When automatic registrationfails, human interaction through a user interface, such as through userinterface function module 230 of FIG. 2, allows for manual registration,through a windows-based interface, such as of the type discussed inconnection with FIG. 6. In this instance a human would manuallyinitially align images more precisely (relying on eyesight and a mousepointer, or, by manually entering via a keyboard the coordinates forpositioning an image) for either pairwise local registration or globalregistration, in order to then allow the computer to automaticallyoptimize the images as discussed herein. The user would not supersedethe system of the present invention, but assist the system to allow forthe automatic convergence to a solution.

[0101] Similarly, in initial calibration, e.g., step 314 in theembodiment of FIG. 3, the number of camera internal and externalparameters is large in the general case. The global optimizationfunction module needs initial conditions in order to converge to theright answer. A user interface, such as of the kind discussed in FIG. 6,provides an interactive tool to initialize those parameters. Thus, inthe event the user has more accurate information about initialparameters than is supplied by the computer default parameters, the usercan input those parameters.

[0102] Human interaction is present throughout to provide feedback tothe computer system of the present invention. The system must have theability to provide feedback in all the nonlinear optimizations to letusers monitor the progress of the system, and allow them to intervenewhen necessary. In a preferred embodiment, the user interface for humaninteraction may be a real-time texture map engine which simulates avirtual camera looking out from the common nodal point of therectilinear images. All images are seen floating in 3D space. The usercan select any image and change its internal and external parametersinteractively in real time. In one preferred embodiment, the images maybe seen as if the user was situated inside of a sphere, termed aspherical coffee-table, illustrated conceptually as virtual realityspace in FIG. 12, with the images 1210 presented tangent to the outsideof the viewing sphere 1220. The arrangement of the images outside thesphere may be arranged by a virtual reality type authoring userinterface, with the user free to arrange which images will be blended,and where, by placing the images at selected tangents on the viewingsphere 1220. Details of this interface can be implemented by one ofordinary skill in the art using the teachings of the present invention.The images, once selected by a user, are pairwise registered,calibrated, globally optimized and blended, as taught herein, toconstruct a panorama that may then be projected onto a sphere, producinga 3D panorama viewing space.

[0103] Although the present invention has been described in terms of thepreferred embodiments above, numerous modifications and/or additions tothe above-described preferred embodiments would be readily apparent toone skilled in the art. It is intended that the scope of the presentinvention extends to all such modifications and/or additions and thatthe scope of the present invention is limited solely by the claims setforth below.

We claim:
 1. A computer system for authoring panoramas, comprising:memory storing data representations of a plurality of images; I/O forinputting and outputting data; function modules in said memory; aprocessor cooperating with said memory and I/O for processinginstructions and data from said memory, said I/O and said functionmodules; a pairwise registration function module in said memory andoperating said processor for pairwise registration of said images, saidpairwise registration function module generating output data related tosaid pairwise registration of said images; a calibration and globaloptimization function module in said memory and operating said processorfor calibration and global optimization of said output data of saidpairwise registration module, said calibration and global optimizationmodule generating output data related to said calibration and globaloptimization of said images; a blending function module in said memoryand operating said processor for generating at least one panorama fromsaid output data from said pairwise registration module and saidcalibration and global optimization module; and a projection functionmodule in said memory and operating said processor for forming apanorama from said images.
 2. The invention according to claim 1,wherein: said pairwise registration function module operates saidprocessor to register a plurality of each of said images that overlapwith one another, by perturbing said images.
 3. The invention accordingto claim 2, wherein: said pairwise registration function module performsa projective transformation between a pair of said images, saidprojective transformation is estimated by said pairwise registrationfunction module by performing a gradient-based optimization thatminimizes a local error function.
 4. The invention according to claim 3,further comprising: a user interface module in said memory, andcooperating with said I/O, operating said processor for initializingsaid projective matrix of said pairwise registration function module,said user interface module operating to process input from a user of thecomputer system, wherein said projective matrix is initialized with saidcamera internal and external parameters from manual input by a user;and, said pairwise registration function projective transformation errorfunction is comprised of the differences in exposures between pixelscomprising the pairs of images, and contains a projective matrixinitialized by camera internal and external parameters relating to saidplurality of images.
 5. The invention according to claim 4, wherein:said internal and external parameters are selected from the groupcomprising: image center position of each of said pairs of images,aspect ratio of said pixels comprising said images, the focal length ofthe camera which took said images; the camera orientation with respectto a common reference frame, camera pan, tilt, roll and skew andbrightness and contrast of said pairs of images.
 6. The inventionaccording to claim 3, wherein: said local error function, eij, is givenby the equation$e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}$

where s_(ij) and b_(ij), the exposure parameters, represent the exposuredifference, I_(i)( ) and I_(j)( ) are pixel intensity values from twoimages, and A_(ij) is the overlap area of said images, the vector X_(i)represents the homogeneous coordinates of the pixel locations making upsaid images, and the matrix M_(ij) represents the matrix that transformscoordinates from image j to image i.
 7. The invention according to claim2, wherein: said pairwise registration function module operates saidprocessor to form a Gaussian pyramid in said memory from said pluralityof overlapping images, said Gaussian pyramid having a plurality oflevels, each subsequent level being finer than each preceding level,said optimization of said local error function is performed onprogressively finer levels of said Gaussian pyramid.
 8. The inventionaccording to claim 7, wherein: the number of pyramid levels of saidGaussian pyramid is computed by said pairwise registration functionmodule according to the formula:$\log_{2}( \frac{l}{l_{\min}} )$

where l is the square root of the smaller eigenvalue of the 2×2 inertialtensor of the overlap polygon region of said plurality of images, andl_(min) is the minimal size of the finest level resolution pyramidlevel.
 9. The invention according to claim 6, wherein: said pairwiseregistration function module operates said processor to form a Gaussianpyramid in said memory from said plurality of overlapping images, saidGaussian pyramid having a plurality of levels, each subsequent levelbeing finer than each preceding level, said optimization of said localerror function is performed on progressively finer levels of saidGaussian pyramid; said pairwise registration function module operatessaid processor to estimate said exposure parameters, s_(ij) and b_(ij),of said local error function by a correlation-based linear search, andattenuates said exposure parameters with progressive damping, saidprogressive damping set to maximum at the coarsest level of said pyramidand decreasing to zero at the finest level.
 10. The invention accordingto claim 2, wherein: a calibration and global optimization functionmodule in said memory and operating said processor for calibration andglobal optimization of said output data of said pairwise registrationmodule, said calibration and global optimization module generatingoutput data related to said calibration and global optimization of allof said overlapping images.
 11. The invention according to claim 10,wherein: said calibration and global optimization function modulecalibrates and optimizes said output data of said pairwise registrationmodule by approximating said local error function.
 12. The inventionaccording to claim 11, wherein said approximation of said local errorfunction, eij, is a Taylor Series expansion given by the formula: e_(ij)(M _(ij))≈e _(ij) ⁰+(M _(ij) −M _(ij) ⁰)^(T) C _(ij)(M _(ij) −M_(ij) ⁰)where e_(ij) ⁰ is a constant representing the minimal valueachieved by said pairwise registration function module; M_(ij) is theprojective matrix of overlapping pairs of images, M_(ij) ⁰ is the vectorrepresenting the optimal projective matrix, and C_(ij) is the Hessianmatrix obtained by the pairwise registration function module whenoptimizing objective function e_(ij).
 13. The invention according toclaim 12, wherein said calibration and global optimization functionmodule calculates a global error function for all overlapping images byapplication of the chain rule to said output data of said pairwiseregistration module, according to the formula:$\frac{\partial E}{\partial( {p_{i},q_{i}} )} = {\sum\limits_{j}{\frac{\partial e_{ij}}{\partial M_{ij}}{\frac{\partial M_{ij}}{\partial( {p_{i},q_{i}} )}.}}}$

and minimizes said global error function, where p_(j)=the vector forcamera internal parameters for each image j_(i) and q_(j)=the cameraorientation with respect to a common reference frame for each image j.14. The invention according to claim 10, further comprising: saidpairwise registration function module operates said processor toestimate a projective transformation for pairwise registration ofimages.
 15. The invention according to claim 14, wherein: saidprojective transformation is initialized by camera internal and externalparameters relating to said plurality of images in said pairwiseregistration function; and said projective transformation is initializedby said output data generated by said calibration and globaloptimization module generating output data, in an iterative manner, and,said calibration and global optimization function module attenuatescamera internal parameters during said calibration and globaloptimization.
 16. The invention according to claim 1, wherein: saidblending function module generating at least one panorama by employingin the image overlap areas of said output data a multi-resolutionweighted average blending.
 17. The invention according to claim 1,further comprising: a user interface function module in said memory andoperating said processor for accepting input from a human user, and foradjusting the areas said images overlap; said blending function modulegenerating at least one panorama by blending said images where theimages overlap, in response to said user interface module.
 18. Theinvention according to claim 1: said blending function module generatessaid panorama from said output data comprising images having imageoverlap regions, said image overlap regions blended by said panoramagenerating function module by performing a grassfire transform on datarelating to said images to generate a blend mask of said image overlapregions.
 19. The invention according to claim 1: said blending functionmodule generates said panorama from said output data comprising imageshaving image overlap regions, said blending function module blending byinstructing said processor to construct Laplacian pyramids stored insaid memory for each plurality of images that overlap in said imageoverlap regions, said Laplacian pyramids having a plurality of levels,said blending function module determining the coarsest level resolutionlevel of the Laplacian pyramids for each of said images, said blendingfunction module generating a blend mask from the image overlap region ofsaid plurality of images, said blending function module constructing aGaussian pyramid from said blend mask for each said plurality of images,and blending said images by multiplying said Laplacian and Gaussianpyramid values for every level of said Laplacian pyramids for every saidplurality of images.
 20. The invention according to claim 1, whereinsaid blending function module blends the high frequency components ofsaid images that overlap over a smaller region in the area of overlapthan the low frequency components of said images.
 21. The inventionaccording to claim 20, wherein said smaller region is on the order of3-5 pixels.
 22. The invention according to claim 1, wherein saidprojection function module forms said panorama on a geometry that isselected from geometries consisting of: cubic, polyhedral, cylindricaland spherical geometries, frusto-conical cones joined at their base withthe apexes pointing away from one another; all geometries employingconical sections, and geometries using the following projections:equidistant, equiangular, ellipsoid, Mercator, transverse Mercator,Oblique Mercator, cylindrical equal-area, Miller cylindrical,equidistant cylindrical, Cassini, Albers equal-area, Lambert conformalconic, equidistant conic, bipolar oblique conic conformal, polyconic,Bonne, azimuthal, orthographic, sterographic, gnomonic, generalperspective, Lambert azimuthal equal-area, azimuthal equidistant,modified-stereographic conformal, space map projections, space obliqueMercator, satellite-tracking projections, pseudocylindrical projections,Van der Grinten, sinusoidal, Mollweide and Eckert IV and VI projections.23. In a computer having a processor which executes instructions storedin memory, a method for generating panoramas comprising the steps of:storing at least two rectilinear images in a computer memory; generatingimage output data by performing pairwise registration of saidrectilinear images as an estimate to a projective transformationrelating to said rectilinear images; generating image output data bycalibration and global optimization of said rectilinear images fromimage output data from said pairwise registration; blending saidrectilinear images to create a panorama from said generated image outputdata from said pairwise registration and said calibration and globaloptimization.
 24. The method according to claim 23, wherein: said stepof blending said rectilinear images to create a panorama comprises thestep of identifying a region of overlap between said rectilinear imagesfrom said generated data; generating a blend mask from said rectilinearimages from said generated data; locating the boundary of said mask byperforming a grassfire transform on said at least two rectilinear imagesindividually.
 25. The method according to claim 23, wherein: said stepof generating data by performing pairwise registration of saidrectilinear images by estimating said projective transformation of saidrectilinear images comprises the steps of: initializing said projectivetransformation with camera internal and external parameters of tworectilinear images that overlap; forming a Gaussian pyramid from saidplurality of overlapping rectilinear images; perturbing said overlappingrectilinear images until a local registration error function isminimized at a first level of said Gaussian pyramid for said images. 26.The method according to claim 25, wherein: wherein said camera internaland external parameters are selected from the group comprising the imagecenter position of each of said pairs of images, aspect ratio of saidpixels comprising said images, the focal length of the camera which tooksaid images; the camera orientation with respect to a common referenceframe, camera pan, tilt, roll and skew and brightness and contrast ofsaid images.
 27. The method according to claim 25, wherein: said localregistration error function, eij, is given by the equation:$e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}$

 where s_(ij) and b_(ij), the exposure parameters, represent theexposure difference, I_(i)( ) and I_(j)( ) are pixel intensity valuesfrom two images, and A_(ij) is the overlap area of said images, thevector X_(i) represents the homogeneous coordinates of the pixellocations making up said images, and the matrix M_(ij) represents thematrix that transforms coordinates from image j to image i; furthercomprising the step of pairwise registration of said overlapping imageson successive iterations of said Gaussian pyramid.
 28. The methodaccording to claim 27, further comprising the step of: attenuating saidestimated exposure parameters, sij and bij, in said equation, at eachsuccessive iteration of said Gaussian pyramid by progressive damping,said damping set to maximum at the coarsest level of said Gaussianpyramid and decreasing to zero at the finest level of said Gaussianpyramid.
 29. The method according to claim 28, further comprising thestep of: estimating said exposure parameters, sij and bij, according toa correlation based linear search.
 30. The method according to claim 27,further comprising the step of: finding said first level of saidGaussian pyramid by successively moving down said Gaussian pyramid untila predetermined number of pixels from each of said images overlap. 31.The method according to claim 30, further comprising the step of:determining the number of levels of said Gaussian pyramid of saidpairwise registration of overlapping images according to the formula:$\log_{2}( \frac{l}{l_{\min}} )$

 where l is the square root of the smaller eigenvalue of the 2×2inertial tensor of the overlap polygon region of said plurality ofimages, and l_(min) is the minimal size of the finest level resolutionpyramid level.
 32. The method according to claim 23, further comprisingthe step of: projecting said panorama formed from said steps of pairwiseregistration and calibration and global optimization, wherein saidprojection function module forms said panorama on a geometry that isselected from geometries consisting of: cubic, polyhedral, cylindricaland spherical geometries, frusto-conical cones joined at their base withthe apexes pointing away from one another; all geometries employingconical sections, and geometries using the following projections:equidistant, equiangular, ellipsoid, Mercator, transverse Mercator,Oblique Mercator, cylindrical equal-area, Miller cylindrical,equidistant cylindrical, Cassini, Albers equal-area, Lambert conformalconic, equidistant conic, bipolar oblique conic conformal, polyconic,Bonne, azimuthal, orthographic, sterographic, gnomonic, generalperspective, Lambert azimuthal equal-area, azimuthal equidistant,modified-stereographic conformal, space map projections, space obliqueMercator, satellite-tracking projections, pseudocylindrical projections,Van der Grinten, sinusoidal, Mollweide and Eckert IV and VI projections.33. The invention according to claim 25, wherein: wherein saidcalibration and global optimization function module calibrates andoptimizes said output data of said pairwise registration module byapproximating said local error function, eij; wherein said step ofestimating said local error function is an approximation of said localerror function given by the formula: e _(ij)(M _(ij))≈e _(ij) ⁰+(M _(ij)−M _(ij) ⁰)^(T) C _(ij)(M _(ij) −M _(ij) ⁰)  where e_(ij) ⁰ is aconstant representing the minimal value achieved by said pairwiseregistration function module; M_(ij) is the projective matrix ofoverlapping pairs of images, M_(ij) ⁰ is the vector representing theoptimal projective matrix, and C_(ij) is the Hessian matrix obtained bythe pairwise registration function module when optimizing objectivefunction e_(ij).
 34. The invention according to claim 33, furthercomprising the steps of: saving said approximated local error functiongenerated during pairwise registration of said images according to thesteps above; computing, for each and every overlapping images, saidcalibration and global optimization of said images, according to theglobal error formula:$\frac{\partial E}{\partial( {p_{i},q_{i}} )} = {\sum\limits_{j}{\frac{\partial e_{ij}}{\partial M_{ij}}\frac{\partial M_{ij}}{\partial( {p_{i},q_{i}} )}}}$

 minimizing said global error function, where p_(j)=the vector forcamera internal parameters for each image j_(i) and q_(j)=the cameraorientation with respect to a common reference frame for each image j.35. The invention according to claim 34, further comprising the stepsof: re-initializing said projective transformation, in an iterativemanner, with said output data generated by said calibration and globaloptimization module.
 36. The method according to claim 23, wherein: saidstep of blending comprises blending high frequency components of saidimages that overlap over a smaller region in the area of overlap thanlow frequency components of said images.
 37. The invention according toclaim 36, wherein said blending for high frequency components of saidimages is performed over said smaller region on the order of 3-5 pixels.38. The invention according to claim 23, wherein said steps furthercomprising: wherein said blending of said rectilinear images to create apanorama comprises the steps of: constructing a Laplacian pyramid forthe overlapping images in an image overlap region where said imagesoverlap, said Laplacian pyramids having a plurality of levels;determining the coarsest level resolution level of the Laplacianpyramids for each of said images; constructing a blend mask from theimage overlap region for said plurality of overlapping images, saidblend mask being a Gaussian pyramid combining said Laplacian pyramid andsaid Gaussian pyramids to blend said images in said image overlapregion.
 39. The method according to claim 38, wherein said step ofconstructing said blend mask is performed by locating the boundary ofsaid blend mask by performing a grassfire transform on at least twooverlapping images.
 40. The invention according to claim 38, furthercomprising the steps of: labeling, during the steps of blending, thepixels comprising each said overlapping images with values to associateeach pixel with a particular image; generating a blend mask from theimage overlap region of said images, and computing a boundary to saidblend mask from said values of said pixel labels and from saidoverlapping images.
 41. The invention according to claim 38, whereinsaid blending function module generating at least one panorama from saidoutput of said pairwise registration module and said calibration andglobal optimization module by employing in the image overlap areas ofsaid output a multi-resolution weighted average blending.
 42. Theinvention according to claim 23, further comprising the steps of:projecting said panorama from said images, after said images have beenpairwise registered, calibrated and globally optimized; projecting saidpanorama onto a view surface geometry that is selected from geometriesconsisting of: cubic, polyhedral, cylindrical and spherical geometries,frusto-conical cones joined at their base with the apexes pointing awayfrom one another; all geometries employing conical sections, andgeometries using the following projections: equidistant, equiangular,ellipsoid, Mercator, transverse Mercator, Oblique Mercator, cylindricalequal-area, Miller cylindrical, equidistant cylindrical, Cassini, Albersequal-area, Lambert conformal conic, equidistant conic, bipolar obliqueconic conformal, polyconic, Bonne azimuthal, orthographic, sterographic,gnomonic, general perspective, Lambert azimuthal equal-area, azimuthalequidistant, modified-stereographic conformal, space map projections,space oblique Mercator, satellite-tracking projections,pseudocylindrical projections, Van der Grinten, sinusoidal, Mollweideand Eckert IV and VI projections.
 43. In a computer having a processorwhich executes instructions stored in memory, a method for generatingpanoramas from two dimensional images, comprising the steps of: storingat least two rectilinear images in a computer memory; blending saidrectilinear images to create a panorama.
 44. The method according toclaim 43, further comprising the steps of: wherein said blending of saidrectilinear images to create a panorama comprises the steps of:constructing a Laplacian pyramid for the overlapping images in an imageoverlap region where said images overlap, said Laplacian pyramids havinga plurality of levels; determining the coarsest level resolution levelof the Laplacian pyramids for each of said images; constructing a blendmask from the image overlap region for said plurality of overlappingimages, said blend mask being a Gaussian pyramid combining saidLaplacian pyramid and said Gaussian pyramids to blend said images insaid image overlap region.
 45. The method according to claim 44, whereinsaid step of constructing said blend mask is performed by locating theboundary of said blend mask by performing a grassfire transform on atleast two overlapping images.
 46. The invention according to claim 44,further comprising the steps of: labeling, during the steps of blending,the pixels comprising each said overlapping images with values toassociate each pixel with a particular image; generating a blend maskfrom the image overlap region of said images, and computing a boundaryto said blend mask from said values of said pixel labels and from saidoverlapping images.
 47. The invention according to claim 44, whereinsaid blending function module generating at least one panorama from saidoutput of said pairwise registration module and said calibration andglobal optimization module by employing in the image overlap areas ofsaid output a multi-resolution weighted average blending.
 48. Theinvention according to claim 43, further comprising the steps of:projecting said panorama from said images, after said images have beenpairwise registered, calibrated and globally optimized; projecting saidpanorama onto a view surface geometry that is selected from geometriesconsisting of: cubic, polyhedral, cylindrical and spherical geometries,frusto-conical cones joined at their base with the apexes pointing awayfrom one another; all geometries employing conical sections, andgeometries using the following projections: equidistant, equiangular,ellipsoid, Mercator, transverse Mercator, Oblique Mercator, cylindricalequal-area, Miller cylindrical, equidistant cylindrical, Cassini, Albersequal-area, Lambert conformal conic, equidistant conic, bipolar obliqueconic conformal, polyconic, Bonne, azimuthal, orthographic,sterographic, gnomonic, general perspective, Lambert azimuthalequal-area, azimuthal equidistant, modified-stereographic conformal,space map projections, space oblique Mercator, satellite-trackingprojections, pseudocylindrical projections, Van der Grinten, sinusoidal,Mollweide and Eckert IV and VI projections.
 49. A computer readablememory storing a program which executes the steps of: storing at least apair of rectilinear images in a computer memory; performing pairwiseregistration of said rectilinear images by calculating a projectivetransformation between said rectilinear images.
 50. The method accordingto claim 49, wherein: said step of performing pairwise registration ofsaid rectilinear images by calculating said projective transformationbetween said rectilinear images comprises the steps of: initializingsaid projective transformation with camera internal and externalparameters of a pair of said rectilinear images that overlap; forming aGaussian pyramid from said pair of overlapping rectilinear images;perturbing said overlapping rectilinear images until a localregistration error function is minimized at a first level of saidGaussian pyramid for said images.
 51. The method according to claim 50,wherein: wherein said camera internal and external parameters areselected from the group comprising the image center position of each ofsaid pairs of images, aspect ratio of said pixels comprising saidimages, the focal length of the camera which took said images; thecamera orientation with respect to a common reference frame, camera pan,tilt, roll and skew and brightness and contrast of said pairs of images.52. The method according to claim 50, wherein: said local registrationerror function, eij, is given by the equation:$e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}$

 where s_(ij) and b_(ij), the exposure parameters, represent theexposure difference, I_(i)( ) and I_(j)( ) are pixel intensity valuesfrom two images, and A_(ij) is the overlap area of said images, thevector X_(i) represents the homogeneous coordinates of the pixellocations making up said images, and the matrix M_(ij) represents thematrix that transforms coordinates from image j to image i; furthercomprising the step of pairwise registration of said overlapping imageson successive iterations of said Gaussian pyramid.
 53. The methodaccording to claim 52, further comprising the step of: attenuating saidestimated exposure parameters, sij and bij, in said equation, at eachsuccessive iteration of said Gaussian pyramid by progressive damping,said damping set to maximum at the coarsest level of said Gaussianpyramid and decreasing to zero at the finest level of said Gaussianpyramid.
 54. The method according to claim 53, further comprising thestep of: estimating said exposure parameters, sij and bij, according toa correlation based linear search.
 55. The method according to claim 52,further comprising the step of: finding said first level of saidGaussian pyramid by iteratively moving down said Gaussian pyramid untila predetermined number of pixels from each of said images overlap. 56.The method according to claim 55, further comprising the step of:determining the number of levels of said Gaussian pyramid of saidpairwise registration of overlapping images according to the formula:${\log \quad}_{2}( \frac{l}{l_{\min}} )$

 where l is the square root of the smaller eigenvalue of the 2×2inertial tensor of the overlap polygon region of said plurality ofimages, and l_(min) is the minimal size of the finest level resolutionpyramid level.
 57. A computer system comprising: a memory storing datarepresentations relating to at least two images; function modules insaid memory; a processor for processing instructions and data from saidmemory and said function modules; a pairwise registration functionmodule in said memory, and operating said processor for pairwiseregistration of said images, said pairwise registration function modulegenerating output data related to said pairwise registration of saidimages.
 58. The invention according to claim 57, wherein: said pairwiseregistration function module operates said processor to register aplurality of each of said images that overlap with one another, byperturbing said images.
 59. The invention according to claim 58,wherein: said pairwise registration function module performs saidpairwise registration by estimating a projective transformation relatinga pair of images that overlap.
 60. The invention according to claim 59,wherein: said projective transformation between said pair of images isestimated by said pairwise registration function module by performing agradient-based optimization that minimizes a local error function thatis comprised of the differences in exposures between pixels comprisingpairs of overlapping images, and contains a projective matrixinitialized by camera internal and external parameters relating to saidplurality of images.
 61. The invention according to claim 60 wherein:said local error function, eij, is given by the equation$e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}$

 where s_(ij) and b_(ij), the exposure parameters, represent theexposure difference, I_(i)( ) and I_(j)( ) are pixel intensity valuesfrom two images, and A_(ij) is the overlap area of said images, thevector X_(i) represents the homogeneous coordinates of the pixellocations making up said images, and the matrix M_(ij) represents thematrix that transforms coordinates from image j to image i, and; saidpairwise registration function module operates said processor to form aGaussian pyramid in said memory from said plurality of overlappingimages, said Gaussian pyramid having a plurality of levels, eachsubsequent level being finer than each preceding level, saidoptimization of said local error function is performed on progressivelyfiner levels of said Gaussian pyramid.
 62. The invention according toclaim 61, wherein: the number of pyramid levels of said Gaussian pyramidis computed by said pairwise registration function module according tothe formula: ${\log \quad}_{2}( \frac{l}{l_{\min}} )$

 where l is the square root of the smaller eigenvalue of the 2×2inertial tensor of the overlap polygon region of said plurality ofimages, and l_(min) is the minimal size of the finest level resolutionpyramid level.
 63. The invention according to claim 61, wherein: saidpairwise registration function module operates said processor toestimate said exposure parameters, s_(ij) and b_(ij) of said local errorfunction by a correlation-based linear search, and attenuates saidexposure parameters with progressive damping, said progressive dampingset to maximum at the coarsest level of said pyramid and decreasing tozero at the finest level.
 64. The invention according to claim 59,further comprising: a user interface module in said memory and operatingsaid processor for initializing said projective transformation of saidpairwise registration function module, said user interface moduleoperating to process input from a user of the computer system, whereinsaid projective transformation is initialized with said camera internaland external parameters from manual input by a user.
 65. The inventionaccording to claim 64, wherein: said internal and external parametersare selected from the group comprising: image center position of each ofsaid pairs of images, aspect ratio of said pixels comprising saidimages, the focal length of the camera which took said images; thecamera orientation with respect to a common reference frame, camera pan,tilt, roll and skew and brightness and contrast of said pairs of images.66. The invention according to claim 64, wherein: said user interfacemodule allows the user the ability to select the arrangement of saidimages.
 67. The invention according to claim 66, wherein: said userinterface module comprises a virtual reality environment comprising aspherical coffee table.
 68. The invention according to claim 67,wherein: said user interface module allows the user to select thearrangement of said images from a plurality of rows of photos on saidspherical coffee table.
 69. The invention according to claim 59,wherein: a calibration and global optimization function module in saidmemory and operating said processor for calibration and globaloptimization of said output data of said pairwise registration module,said calibration and global optimization module generating output datarelated to said calibration and global optimization of all of saidoverlapping images.
 70. The invention according to claim 69, wherein:wherein said calibration and global optimization function modulecalibrates and optimizes said output data of said pairwise registrationmodule by approximating a local error function.
 71. The inventionaccording to claim 70, wherein said approximation of said local errorfunction is a Taylor Series expansion given by the formula: e _(ij)(M_(ij))≈e _(ij) ⁰+(M _(ij) −M _(ij) ⁰)^(T) C _(ij)(M _(ij) −M _(ij)⁰)where e_(ij) ⁰ is a constant representing the minimal value achievedby said pairwise registration function module; M_(ij) is the projectivematrix of overlapping pairs of images, M_(ij) ⁰ is the vectorrepresenting the optimal projective matrix, and C_(ij) is the Hessianmatrix obtained by the pairwise registration function module whenoptimizing objective function e_(ij).
 72. The invention according toclaim 70, wherein said calibration and global optimization functionmodule calculates a global error function for all overlapping images byapplication of the chain rule to said output data of said pairwiseregistration module, according to the formula:$\frac{\partial E}{\partial( {p_{i},q_{i}} )} = {\sum\limits_{j}{\frac{\partial e_{ij}}{\partial M_{ij}}\frac{\partial M_{ij}}{\partial( {p_{i},q_{i}} )}}}$

and minimizes said global error function, where p_(j)=the vector forcamera internal parameters for each image j, and q_(j)=the cameraorientation with respect to a common reference frame for each image j.73. The invention according to claim 69, wherein: said projectivetransformation initialized by camera internal and external parametersrelating to said plurality of images in said pairwise registrationfunction is initialized by said output data generated by saidcalibration and global optimization module generating output data, in aniterative manner.
 74. The invention according to claim 73, wherein: saidcalibration and global optimization function module attenuates camerainternal parameters during said calibration and global optimization. 75.The invention according to claim 69, further comprising: a blendingfunction module in said memory and operating said processor forgenerating at least one panorama from said output data from saidpairwise registration module and said calibration and globaloptimization module.
 76. The invention according to claim 75, wherein:said blending function module generating at least one panorama from saidoutput of said pairwise registration module and said calibration andglobal optimization module by employing in the image overlap areas ofsaid output a multi-resolution weighted average blending.
 77. Theinvention according to claim 76, further comprising: a user interfacemodule in said memory and operating said processor for accepting inputfrom a human user, and for selecting said images and adjusting the areassaid images overlap; said blending function module generating at leastone panorama from said output of said pairwise registration module andsaid calibration and global optimization module by employing in theimage overlap areas of said output a weighted average blending, inresponse to said user interface module.
 78. The invention according toclaim 76, further comprising: a user interface module in said memory andoperating said processor for accepting input from a human user, and forselecting said images and adjusting the areas. said images overlap,wherein said user interface module allows the user the ability to selectthe arrangement of said images.
 79. The invention according to claim 78,wherein: said user interface module comprises a virtual realityenvironment comprising a spherical coffee table.
 80. The inventionaccording to claim 79, wherein: said user interface module allows theuser to select the arrangement of said images from a plurality of rowsof photos.
 81. The invention according to claim 76, wherein: saidblending function module generates said panorama from said output datacomprising images having image overlap regions, said image overlapregions blended by said panorama generating function module byperforming a grassfire transform on data relating to said images togenerate a blend mask of said images.
 82. The invention according toclaim 76, wherein said blending function module blends high frequencycomponents of said images that overlap over a smaller region in the areaof overlap than low frequency components of said images.
 83. Theinvention according to claim 82, wherein said smaller region is on theorder of 3-5 pixels.
 84. The invention according to claim 76, wherein:said blending function module generates said panorama from said outputdata comprising images having image overlap regions, said blendingfunction module blending by instructing said processor to constructLaplacian pyramids stored in said memory for each plurality of imagesthat overlap in said image overlap regions, said Laplacian pyramidshaving a plurality of levels, said blending function module determiningthe coarsest level resolution level of the Laplacian pyramids for eachof said images, said blending function module generating a blend maskfrom the image overlap region of said plurality of images, said blendingfunction module constructing a Gaussian pyramid from said blend mask foreach said plurality of images, and blending said images by multiplyingsaid Laplacian and Gaussian pyramid values for every level of saidLaplacian pyramids for every said plurality of images.
 85. The inventionaccording to claim 76, wherein: said blending function labels the pixelscomprising each said overlapping images with values to associate eachpixel with a particular image. said blending function module generatinga blend mask from the image overlap region of said images, and computinga boundary to said blend mask from said values of said pixel labels andfrom said overlapping images.
 86. The invention according to claim 76,further comprising: a projection function module in said memory andoperating said processor for forming a panorama from said images, saidprojection function forms said panorama on a geometry that is selectedfrom geometries consisting of: cubic, polyhedral, cylindrical andspherical geometries, frusto-conical cones joined at their base with theapexes pointing away from one another; all geometries employing conicalsections, and geometries using the following projections: equidistant,equiangular, ellipsoid, Mercator, transverse Mercator, Oblique Mercator,cylindrical equal-area, Miller cylindrical, equidistant cylindrical,Cassini, Albers equal-area, Lambert conformal conic, equidistant conic,bipolar oblique conic conformal, polyconic, Bonne, azimuthal,orthographic, sterographic, gnomonic, general perspective, Lambertazimuthal equal-area, azimuthal equidistant, modified-stereographicconformal, space map projections, space oblique Mercator,satellite-tracking projections, pseudocylindrical projections, Van derGrinten, sinusoidal, Mollweide and Eckert IV and VI projections.
 87. Acomputer system comprising: a memory storing data representationsrelating to at least a pair of overlapping images; a blending functionmodule in said memory; a processor for processing instructions and datafrom said memory and said blending function module; said blendingfunction module in said memory operating said processor for blendingsaid images to form a panorama from said images.
 88. The inventionaccording to claim 87, further comprising: said blending function modulegenerating at least one panorama by employing in the image overlap areasof said images a multi-resolution weighted average blending.
 89. Theinvention according to claim 87, further comprising: a user interfacemodule in said memory and operating said processor for accepting inputfrom a human -user, and for selecting said images and adjusting theareas said images overlap; said blending function module generating atleast one panorama by blending said image overlap areas, in response tosaid user interface module.
 90. The invention according to claim 87,wherein: said blending function module generates said panorama byperforming a grassfire transform on data relating to said images togenerate a blend mask of said images.
 91. The invention according toclaim 87, wherein: said blending function module generates said panoramaby instructing said processor to construct Laplacian pyramids stored insaid memory for each plurality of images that overlap in image overlapregions, said Laplacian pyramids having a plurality of levels, saidblending function module determining the coarsest level resolution levelof the Laplacian pyramids for each of said images, said blendingfunction module generating a blend mask from the image overlap region ofsaid plurality of images, said blending function module constructing aGaussian pyramid from said blend mask for each said plurality of images,and blending said images by multiplying said Laplacian and Gaussianpyramid values for every level of said Laplacian pyramids for every saidplurality of images.
 92. The invention according to claim 87, furthercomprising: a projection function module in said memory and operatingsaid processor for forming a panorama from said images, said projectionfunction forms said panorama on a geometry that is selected fromgeometries consisting of: cubic, polyhedral, cylindrical and sphericalgeometries, frusto-conical cones joined at their base with the apexespointing away from one another; all geometries employing conicalsections, and geometries using the following projections: equidistant,equiangular, ellipsoid, Mercator, transverse Mercator, Oblique Mercator,cylindrical equal-area, Miller cylindrical, equidistant cylindrical,Cassini, Albers equal-area, Lambert conformal conic, equidistant conic,bipolar oblique conic conformal, polyconic, Bonne, azimuthal,orthographic, sterographic, gnomonic, general perspective, Lambertazimuthal equal-area, azimuthal equidistant, modified-stereographicconformal, space map projections, space oblique Mercator,satellite-tracking projections, pseudocylindrical projections, Van derGrinten, sinusoidal, Mollweide and Eckert IV and VI projections.
 93. Theinvention according to claim 87, further comprising: a user interfacemodule in said memory and operating said processor, said user interfacemodule operating to process input from a user of the computer system,said user interface module allowing the user the ability to manuallyselect the arrangement of said images to be blended.
 94. The inventionaccording to claim 93, wherein: said user interface module comprises avirtual reality environment comprising a spherical coffee table.
 95. Theinvention according to claim 94, wherein: said user interface moduleallows the user to select the arrangement of said images from aplurality of rows or columns of photos.
 96. The invention according toclaim 87, wherein: said blending function labels the pixels comprisingeach said overlapping images with values to associate each pixel with aparticular image, said blending function module generating a blend maskfrom the image overlap region of said. images, and computing a boundaryto said blend mask from said values of said pixel labels and from saidoverlapping images.
 97. The invention according to 96, wherein saidblend mask is constructed by said blending module function performing agrassfire transform on said overlapping images.
 98. The inventionaccording to claim 87, further comprising: a pairwise registrationfunction module in said memory, and operating said processor forpairwise registration of said images, said pairwise registrationfunction module generating output data related to said pairwiseregistration of said images, said pairwise registration function moduleoperates said processor to register a plurality of each of said imagesthat overlap with one another, by perturbing said images.
 99. Theinvention according to claim 98, wherein: said pairwise registrationfunction module performs a pairwise registration by performing aprojective transformation between said images, said projectivetransformation is estimated by said pairwise registration functionmodule by performing a gradient-based optimization that minimizes alocal error function.
 100. The invention according to claim 99, furthercomprising: a user interface module in said memory and operating saidprocessor for initializing said projective matrix of said pairwiseregistration function module, said user interface module operating toprocess input from a user of the computer system, wherein saidprojective matrix is initialized with said camera internal and externalparameters from manual input by a user; and, said error function of saidprojective transformation is comprised of the differences in exposuresbetween pixels comprising said pairs of overlapping images, and containsa projective matrix initialized by camera internal and externalparameters relating to said plurality of images.
 101. The inventionaccording to claim 100, wherein: said internal and external parametersare selected from the group comprising: image center position of each ofsaid pairs of images, aspect ratio of said pixels comprising saidimages, the focal length of the camera which took said images; thecamera orientation with respect to a common reference frame, camera pan,tilt, roll and skew and brightness and contrast of said pairs of images.102. The invention according to claim 100, wherein: said local errorfunction, eij, is given by the equation$e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}$

 where s_(ij) and b_(ij), the exposure parameters, represent theexposure difference, I_(i)( ) and I_(j)( ) are pixel intensity valuesfrom two images, and A_(ij) is the overlap area of said images, thevector X_(i) represents the homogeneous coordinates of the pixellocations making up said images, and the matrix M_(ij) represents thematrix that transforms coordinates from image j to image i.
 103. Theinvention according to claim 102, wherein: said pairwise registrationfunction module operates said processor to form a Gaussian pyramid insaid memory from said plurality of overlapping images, said Gaussianpyramid having a plurality of levels, each subsequent level being finerthan each preceding level, said optimization of said local errorfunction is performed on progressively finer levels of said Gaussianpyramid; and said pairwise registration function module operates saidprocessor to estimate said exposure parameters, s_(ij) and b_(ij), ofsaid local error function by a correlation-based linear search, andattenuates said exposure parameters with progressive damping, saidprogressive damping set to maximum at the coarsest level of said pyramidand decreasing to zero at the finest level.
 104. The invention accordingto claim 98, wherein: said pairwise registration-function moduleoperates said processor to form a Gaussian pyramid in said memory fromsaid plurality of overlapping images, said Gaussian pyramid having aplurality of levels, each subsequent level being finer than eachpreceding level, said optimization of said local error function isperformed on progressively finer levels of said Gaussian pyramid. 105.The invention according to claim 104, wherein: the number of pyramidlevels of said Gaussian pyramid is computed by said pairwiseregistration function module according to the formula:${\log \quad}_{2}( \frac{l}{l_{\min}} )$

10 where l is the square root of the smaller eigenvalue of the 2×2inertial tensor of the overlap polygon region of said plurality ofimages, and l_(min) is the minimal size of the finest level resolutionpyramid level.
 106. The invention according to claim 98, wherein: acalibration and global optimization function module in said memory andoperating said processor for calibration and global optimization of saidoutput data of said pairwise registration module, said calibration andglobal optimization module generating output data related to saidcalibration and global optimization of all of said overlapping images.107. The invention according to claim 106, wherein: wherein saidcalibration and global optimization function module calibrates andoptimizes said output data of said pairwise registration module byapproximating a local error function.
 108. The invention according toclaim 107, wherein said approximation of said local error function is aTaylor Series expansion given by the formula: e _(ij)(M _(ij))≈e _(ij)⁰+(M _(ij) −M _(ij) ⁰)^(T) C _(ij)(M _(ij) −M _(ij) ⁰)where e_(ij) ⁰ isa constant representing the minimal value achieved by said pairwiseregistration function module; M_(ij) is the projective matrix ofoverlapping pairs of images, M_(ij) ⁰ is the vector representing theoptimal projective matrix, and C_(ij) is the Hessian matrix obtained bythe pairwise registration function module when optimizing objectivefunction e_(ij).
 109. The invention according to claim 108, wherein saidcalibration and global optimization function module calculates a globalerror function for all overlapping images by application of the chainrule to said output data of said pairwise registration module, accordingto the formula:$\frac{\partial E}{\partial( {p_{i},q_{i}} )} = {\sum\limits_{j}{\frac{\partial e_{ij}}{\partial M_{ij}}\frac{\partial M_{ij}}{\partial( {p_{i},q_{i}} )}}}$

and minimizes said global error function, where p_(j)=the vector forcamera internal parameters for each image j, and q_(j)=the cameraorientation with respect to a common reference frame for each image j.110. The invention according to claim 106, wherein: said pairwiseregistration function module performs a pairwise registration byperforming a projective transformation between said images, saidprojective transformation- is estimated by said pairwise registrationfunction module by performing a gradient-based optimization thatminimizes a local error function; and the projective transformationbetween said pair of said images is initialized by camera internal andexternal parameters relating to said images in said pairwiseregistration function module by said calibration and global optimizationmodule, in an iterative manner.
 111. The invention according to claim110, wherein: said calibration and global optimization function moduleattenuates camera internal parameters during said calibration and globaloptimization.
 112. The invention according to claim 106, furthercomprising: a blending function module in said memory and operating saidprocessor for generating at least one panorama from said output datafrom said pairwise registration module and said calibration and globaloptimization module.
 113. An apparatus for generating a panoramacomprising: means for pairwise registration of images having overlappingregions; means for calibration and global optimization of said images;means for blending said images after pairwise registration, calibrationand global optimization of said images; means for projecting saidblended images into a panorama.
 114. The invention according to claim113, wherein: said pairwise registration means operates said processorto register a plurality of each of said images that overlap with oneanother, through perturbation of said images.
 115. The inventionaccording to claim 114, wherein: said pairwise registration meansestimates a projective transformation between a pair of images byperforming a gradient-based optimization that minimizes a local errorfunction.
 116. The invention according to claim 115, further comprising:a user interface module in said memory and operating said processor forinitializing said projective matrix of said pairwise registration means,said user interface module operating to process input from a user of thecomputer system, wherein said projective matrix is initialized with saidcamera internal and external parameters from manual input by a user; andthe local error function of said pairwise registration means iscomprised of the differences in exposures between pixels comprising thepairs of images.
 117. The invention according to claim 116, wherein:said projective transformation contains a projective matrix parametrizedby said internal and external parameters, selected from the groupcomprising: image center position of each of said pairs of images,aspect ratio of said pixels comprising said images, the focal length ofthe camera which took said images; the camera orientation with respectto a common reference frame, camera pan, tilt, roll and skew andbrightness and contrast of said pairs of images, and, said projectivetransformation is initialized by camera internal and external parametersrelating to the plurality of images.
 118. The invention according toclaim 115, wherein: said local error function, eij, is given by theequation$e_{ij} = {\frac{1}{A_{ij}}{\sum\limits_{overlap}( {{s_{ij}{I_{j}( X_{j} )}} + b_{ij} - {I_{i}( {M_{ij}X_{j}} )}} )^{2}}}$

 where s_(ij) and b_(ij), the exposure parameters, represent theexposure difference, I_(i)( ) and I_(j)( ) are pixel intensity valuesfrom two images, and A_(ij) is the overlap area of said images, thevector X_(i) represents the homogeneous coordinates of the pixellocations making up said images, and the matrix M_(ij) represents thematrix that transforms coordinates from image j to image i.
 119. Theinvention according to claim 118, wherein: said pairwise registrationmeans operates said processor to form a Gaussian pyramid in said memoryfrom said plurality of overlapping images, said Gaussian pyramid havinga plurality of levels, each subsequent level being finer than eachpreceding level, said optimization of said local error function isperformed on progressively finer levels of said Gaussian pyramid. 120.The invention according to claim 119, wherein: the number of pyramidlevels of said Gaussian pyramid is computed by said pairwiseregistration means according to the formula:$\log_{2}( \frac{l}{l_{\min}} )$

 where l is the square root of the smaller eigenvalue of the 2×2inertial tensor of the overlap polygon region of said plurality ofimages, and l_(min) is the minimal size of the finest level resolutionpyramid level.
 121. The invention according to claim 119, wherein: saidpairwise registration means operates said processor to estimate saidexposure parameters, s_(ij) and b_(ij), of said local error function bya correlation-based linear search, and attenuates said exposureparameters with progressive damping, said progressive damping set tomaximum at the coarsest level of said pyramid and decreasing to zero atthe finest level.
 122. The invention according to claim 115, wherein: acalibration and global optimization means in said memory and operatingsaid processor for calibration and global optimization of said outputdata of said pairwise registration module, said calibration and globaloptimization module generating output data related to said calibrationand global optimization of all of said overlapping images.
 123. Theinvention according to claim 122, wherein: wherein said calibration andglobal optimization means calibrates and optimizes said output data ofsaid pairwise registration module by approximating said local errorfunction, eij.
 124. The invention according to claim 123, wherein saidapproximation of said local error function is a Taylor Series expansiongiven by the formula: e _(ij)(M _(ij))≈e _(ij) ⁰+(M _(ij) −M _(ij)⁰)^(T) C _(ij)(M _(ij) −M _(ij) ⁰)where e_(ij) ⁰ is a constantrepresenting the minimal value achieved by said pairwise registrationmeans; M_(ij) is the projective matrix of overlapping pairs of images,M_(ij) ⁰ is the vector representing the optimal projective matrix, andC_(ij) is the Hessian matrix obtained by the pairwise registration meanswhen optimizing objective function e_(ij).
 125. The invention accordingto claim 123, wherein said calibration and global optimization meanscalculates a global error function for all overlapping images byapplication of the chain rule to said output data of said pairwiseregistration module, according to the formula:$\frac{\partial E}{\partial( {p_{i},q_{i}} )} = {\sum\limits_{j}{\frac{\partial e_{ij}}{\partial M_{ij}}\frac{\partial M_{ij}}{\partial( {p_{i},q_{i}} )}}}$

and minimizes said global error function, where p_(j)=the vector forcamera internal parameters for each image j, and q_(j)=the cameraorientation with respect to a common reference frame for each image j.126. The invention according to claim 122, wherein: said projectivetransformation is initialized by camera internal and external parametersrelating to said plurality of images in said pairwise registrationfunction, said initialization is by said output data generated by saidcalibration and global optimization module, and said initialization isperformed by said calibration and global optimization module in aniterative manner.
 127. The invention according to claim 126, wherein:said calibration and global optimization means attenuates camerainternal parameters during said calibration and global optimization.128. The invention according to claim 122, further comprising: ablending means in said memory and operating said processor forgenerating at least one panorama from said output data from saidpairwise registration module and said calibration and globaloptimization module.
 129. The invention according to claim 128, wherein:said blending means generating at least one panorama from said output ofsaid pairwise registration module and said calibration and globaloptimization module by employing in the image overlap areas of saidoutput a multi-resolution weighted average blending.
 130. The inventionaccording to claim 129, further comprising: a user interface module insaid memory and operating said processor for accepting input from ahuman user, and for selecting said images and adjusting the areas saidimages overlap.
 131. The invention according to claim 129, furthercomprising: a projecting means in said memory and operating saidprocessor for forming a panorama from said images, said projectionfunction forms said panorama on a geometry that is selectedfrom-geometries consisting of: cubic, polyhedral, cylindrical andspherical geometries, frusto-conical cones joined at their base with theapexes pointing away from one another; all geometries employing conicalsections, and geometries using the following projections: equidistant,equiangular, ellipsoid, Mercator, transverse Mercator, Oblique Mercator,cylindrical equal-area, Miller cylindrical, equidistant cylindrical,Cassini, Albers equal-area, Lambert conformal conic, equidistant conic,bipolar oblique conic conformal, polyconic, Bonne, azimuthal,orthographic, sterographic, gnomonic, general perspective, Lambertazimuthal equal-area, azimuthal equidistant, modified-stereographicconformal, space map projections, space oblique Mercator,satellite-tracking projections, pseudocylindrical projections, Van derGrinten, sinusoidal, Mollweide and Eckert IV and VI projections. 132.The invention according to claim 129, wherein: said blending meansgenerates said panorama from said output data comprising images havingimage overlap regions, said image overlap regions blended by saidpanorama generating means by performing a grassfire transform on datarelating to said images to generate a mask of said images.
 133. Theinvention according to claim 129, wherein: said blending means generatessaid panorama from said output data comprising images having imageoverlap regions, said blending means blending by instructing said-processor to construct Laplacian pyramids stored in said memory for eachplurality of images that overlap in said image overlap regions, saidLaplacian pyramids having a plurality of levels, said blending meansdetermining the coarsest level resolution level of the Laplacianpyramids for each of said images, said blending means generating a blendmask from the image overlap region of said plurality of images, saidblending means constructing a Gaussian pyramid from said blend mask foreach said plurality of images, and blending said images by multiplyingsaid Laplacian and Gaussian pyramid values for every level of saidLaplacian pyramids for every said plurality of images.
 134. Theinvention according to claim 129, wherein: said blending function labelsthe pixels comprising each said overlapping images with values toassociate each pixel with a particular image, said blending meansgenerating a blend mask from the image overlap region of said images,and computing a boundary to said blend mask from said values of saidpixel labels and from said overlapping images.
 135. A computer readablememory comprising: a memory storing a first data structure for pairwiseprojective registration of overlapping images; a memory storing a seconddata structure for calibration and global optimization of data from saidfirst data structure; a memory storing a third data structure forblending data from said first and second data structures; whereby apanorama of three-dimensional images may be constructed fromtwo-dimensional images from said data in said first, second and thirddata structures.
 136. The invention according to claim 135, wherein:said first data structure comprises parameters in a projectiveregistration matrix; said second data structure comprises a constant,e_(ij) representing the minimal value achieved in said projectiveregistration of overlapping images, a vector representing the optimalprojective matrix obtained for said constant e_(ij), and a Hessianmatrix.